Precision modeling of redshift-space distortions from multi-point propagator 

expansion 
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Using a full implementation of resummed perturbation theory (PT) from a multi-point propagator 
expansion, we put forward new theoretical predictions for the two-point statistics of matter fluctua- 
tions in redshift space. The predictions consistently include PT corrections up to the two-loop order 
and are based on an improved prescription of the redshift-space distortions (RSD) that properly 
takes into account their non-Gaussian impact from a systematic low-fc expansion. In contrast to the 
previous studies that partly used standard PT calculations, the present treatment is able to provide 
a consistent prediction for both power spectra and correlation functions. These results are compared 
with TV-body simulations with which a very good agreement is found up to the quadrupole moment. 
The theoretical predictions for the hexadecapole moment of the power spectra are however found 
to significantly departs from the numerical results at low redshift. We examine this issue and found 
it to be likely related to an improper modeling of the RSD damping effects on which this moment 
shows large dependence. 

PACS numbers: 98.80.-k 



I. INTRODUCTION 



The large-scale structure of the Universe observed via 
galaxy redshift surveys appears distorted due to the pe- 
culiar velocity of galaxies, known as the redshift-space 
distortions (RSD) (e.g., [10). The RSD breaks the sta- 
tistical isotropy, and the resultant galaxy clustering ex- 
hibits characteristic anisotropics along observer's line-of- 
sight direction by the two competitive effects, i.e., Kaiser 
and Finger-of-God effects While the latter effect 

mostly comes from the virialized random motion of the 
mass (or galaxy) in halos, the Kaiser effect apparently 
enhances the clustering amplitude along a line-of-sight 
direction, and the effect is dominated by the large-scale 
coherent motion. In particular, the strength of the Kaiser 
effect is simply described by the linear theory, and is 
characterized by the growth-rate parameter, /, defined 
as / = d In D + /d In a, where the quantities D + and a are 
the linear growth factor and scale factor of the Universe, 
respectively (e.g., 0,@])- Thus, the Kaiser effect can be 
used as a useful tool to measure the growth of cosmologi- 
cal structure, and combining the distance measurements, 
the measurement of RSD offers a unique opportunity to 
test the gravity theory on cosmological scales (e.g., 0- 
[l3j). Note that the galaxy redshift surveys also provide 
a way to measure the cosmological distance, utilizing the 
so-called Alcock-Paczynski effect [l4j]. With the baryon 
acoustic oscillations (BAOs) imprinted on the large-scale 
structure as a robust standard ruler, we can thus measure 
the angular diameter distance Da{z) and the Hubble pa- 
rameter H(z) of distant galaxies at redshift z separately 
through the Alcock-Paczynski effect (e.g., |15l - [l9| ). 



Hence, the large-scale galaxy distribution offers a dual 
cosmological probe to simultaneously constrain both the 
cosmic expansion and structure growth, from which we 
can directly test gravity, and even address the nature 
of late-time cosmic acceleration (see [H I20l - l23j for re- 
cent status). This is the main reason why ongoing 
and planned galaxy surveys aim at precisely measuring 
the RSD and BAOs through the clustering statistics of 
galaxy distribution. On the other hand, a pursuit of such 
a precision measurement poses several challenging issues 
in theory of large-scale structure. One important issue 
is the development of precision theoretical tools to com- 
pute the clustering statistics of the large-scale structure. 
While the RSD and BAO are measured from the galaxy 
clustering at the scales close to the linear regime of the 
gravitational evolution, nonlinearity of both the gravity 
and the RSD is known to play a crucial role in precise 
estimate of the parameters /, Da, and H (e.g., [2#l2q ). 

The aim of the present paper is to address such an issue 
and to investigate the extent to which we can accurately 
compute the redshift-space power spectrum and correla- 
tion function based on the perturbation theory (PT). In 
redshift space, a key point is that the applicable range 
of linear theory prediction is fairly narrower than that in 
real space, and the corrections coming from the nonlin- 
ear clustering and RSD need to be properly incorporated 
into theoretical predictions even on such large-scales as 
k < 0.1 /iMpc -1 . This is because, in terms of real-space 
quantities, the redshift-space power spectrum and/or cor- 
relation function cannot simply be expressed as the large- 
scale two-point statistics of the underlying fields and are 
significantly affected by the small scale physics@- Thus, 
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for a reliable theoretical predictions with a wider appli- 
cable range, a sophisticated treatment of both the grav- 
itational clustering and RSD is required. 

In our previous studies, we have proposed an im- 
proved model of RSD relevant in the weakly nonlinear 
regime [24], HH (see [27l - [3lT | for other improved models or 
treatments). The model properly accounts of the non- 
Gaussian nature of the RSD based on the low-A: expan- 
sion. The resulting power spectrum expression differs 
from the one of the so-called streaming model frequently 
used in the literature (e.g., [2j, |32], |33j) with the in- 
troduction of additional corrections associated with non- 
linear couplings between velocity and density fields [2~3 |. 
A detailed investigation revealed that these corrections 
can give an important contribution to the acoustic struc- 
ture of BAOs which gives rise to a slight increase in the 
amplitude of monopole and quadrupole spectra. While 
the model has been originally proposed for the matter 
power spectrum, with an improved PT of gravitational 
clustering as well as an appropriate parameterization of 
galaxy/halo bias, it has been shown to successfully de- 
scribe not only the matter but also the halo power spectra 
in TV-body simulations [261 ] . 

In this paper, the third of a series on this method, we 
present consistent PT calculations aiming at giving pre- 
dictions in both Fourier and configuration spaces. In all 
of our previous works, the PT predictions have been done 
in some heterogeneous ways. That is, while the improved 
PT scheme by Ref. [H, HH is applied to compute the 
power spectra for density and velocity fields, the stan- 
dard PT treatment is partly used to compute the new 
corrections of the RSD, for which we have only given the 
predictions based on the leading-order PT calculation. It 
is generally known that the standard PT treatment pro- 
duces an ill-behaved PT expansion that has poor con- 
vergence properties (e.g., |35l - l37| ). This is indeed the 
case when we compute the redshift-space power spec- 
trum. Because of the bad UV behavior, the correction 
terms computed with standard PT eventually become 
dominant at small scales. Then, a naive computation of 
the correlation function through the direct integration of 
the power spectrum ceases to give a physically reason- 
able answer. This is one of the main reasons why the 
previous works have focused on the redshift-space power 
spectrum, not the correlation function. 

Here, in order to remedy the bad UV behaviors, we 
will apply the specific resummcd PT scheme referred to 
as the multi-point propagator expansion or T expansion 
[38| . The advantage of this scheme is that the non- 
perturbative properties, which can be obtained in stan- 
dard PT by summing up infinite series of PT expansions, 
are whole encapsulated in the multi-point propagators, 
with which any statistical quantities including the power 
spectrum, bispectrum and trispectrum can be built up. 
In particular, the asymptotic behaviors of the multi-point 
propagators are analytically known [38l . [39j . and with 
a help of these analytic properties, a novel regularized 
treatment has been proposed [40|, allowing us to consis- 



tently interpolate between standard PT results at low-fc 
and the expected resummed behavior at high-fc. In this 
paper, making full use of the regularized V expansion, 
we are now able to give a consistent calculation for both 
the power spectrum and correlation function in redshift 
space. We will present the results including the PT cor- 
rections up to the next-to-next-to-leading order, i.e. two- 
loop order, and compare those with TV-body simulations 
in detail. With the calculations at the two-loop order, 
we will also discuss the potential impact of the higher- 
order non-Gaussian corrections coming from the RSD. 
While the model of RSD has been derived based on the 
low-fc expansion, one of the corrections, comparable to 
the two-loop order in PT expansion, has been ignored in 
previous studies without any justification. Here, we will 
explicitly quantify the size of this correction, and with 
a help of the iV-body simulations, we will explore the 
potential systematics of our model predictions. 

The paper is organized as follows. In Sec. UH we begin 
by briefly reviewing the model of RSD. Sec. IIIII describes 
the multi-point propagator expansion which we will apply 
to the predictions of redshift-space power spectrum and 
correlation function. Based on the regularized scheme, 
the basic formalism to compute the propagators and the 
statistical quantities built up with those is presented, to- 
gether with technical detail of the implementation. Then, 
Sec.|lV]presents the main results. The PT predictions up 
to the two-loop order are compared with iV-body simula- 
tions. In Sec. El the impact of higher-order corrections of 
RSD is discussed, and the correction that has been so far 
neglected is explicitly computed. With a help of ./V-body 
simulations, a potential systematics in the model predic- 
tion is also investigated. Finally, Sec. [VI] is devoted to 
the summary and conclusions. 



II. THE MODEL OF REDSHIFT-SPACE 
DISTORTIONS 

Here, we briefly review the model prescription of 
redshift-space distortions (RSD) proposed by Ref. [24j . 
and discuss several remarks on the perturbation theory 
(PT) treatment in computing redshift-space power spec- 
trum or correlation function. 

We begin by writing the exact expression for redshift- 
space power spectrum. Let us denote the density and 
velocity fields by 5 and v. Owing to the distant-observer 
approximation, which is usually valid for the observation 
of distant galaxies of our interest, one can write (e.g., 

mm mi) 

P( s \k)= J d 3 xe ik - x (e- ik >*f Au » 

x {S(r) + fW z u z (r)} {S(r') + fV z u z (r')}), (1) 

where x = r — r' denotes the separation in real space 
and (■ ■ • ) indicates an ensemble average. In the above 
expression, the z-axis is taken as observer's line-of-sight 
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direction, and we define the directional cosine [i by // = 
k z /k. Further, we defined u z (r) = —v z (r)/(aHf), and 
Au z — u z (r) — u z (r') for the line-of-sight component of 
the velocity field. Here, the function / is the logarithmic 
derivative of the linear growth function D(z) given by 
/ = dlxiD(z)/dlxia. Note that the above expression has 
been derived without invoking the dynamical information 
for velocity and density fields, i.e., the Euler equation 
and/or continuity equation. 

Clearly, the redshift-space power spectrum p( s \ given 
as the function of wavenumber k and the directional co- 
sine fj,, cannot be derived from the mere two-point statis- 
tics of the underlying fields. If we expand the exponent 
in the bracket, it can be described by the infinite series 
of multi-point spectra of density and velocity fields. In 
order to derive an analytically tractable expression ap- 
plicable to the large-scale cosmological observations, we 
rewrite the expression (TTJ) with 



P (S W) = Jd 3 xe ik - X (^A 2 A 3 ), 



(2) 



where the quantities j\, Ai (i — 1,2,3) are respectively 
given by 

ji=-i kfif, 

A! = u z (r) - u z (r'), 

A 2 =5{r) + fV z u z (r), 

A 3 =5(r')+ fV z u z (r'). (3) 

We shall rewrite the ensemble average (e-? 1 1 A 2 A 3 ) in 
terms of cumulants. To do this, we use the relation be- 
tween the cumulant and moment generating functions. 
For the stochastic vector field A = {Ai,A2, A3}, we have 
(e.g., 0,13): 



(e3' A ) = exp 



(4) 



with j being an arbitrary constant vector, j — 
{31,32,33}- Taking the derivative twice with respect to 
32 and j'3, we then set 32 = 33 = 0. We obtain @ 

(e^ M A 2 A 3 ) =exp{(e^) c } 

x [(e^A 2 A 3 ) c + (e^A 2 ) c (e^A 3 ) c ] . (5) 

Substituting this into Eq.@, we arrive at 



P< s >(fc,M) = j d a x e ikx exp{(e^i} c } 

x [(e^A 2 A 3 ) c + (e^A 2 ) c (e^A 3 ) c 



(6) 



This expression clearly reveals the coupling between 
density and velocity fields associated with the Kaiser 
and Finger-of-God effects. In addition to the prefactor 
exp {(e^ lAl )c}> the ensemble averages over the quanti- 
ties A 2 and A 3 responsible for the Kaiser effect all in- 
clude the exponential factor eP 1 1 , which can produce a 
non-negligible correlation between density and velocity. 

Based on Eq. © , the proposition raised by Ref. [24| is 
the following: 



• The exponential prefactor, exp {{ei lAl ) c } , which is 
responsible for the suppression of power spectrum 
due to the virialized random and coherent motion, 
turns out to mainly affect the broadband shape of 
the power spectrum (Sec.III-B of [24j). Neverthe- 
less, the effect of this is known to be significant on 
BAO scales, and seems difficult to treat perturba- 
tively. Thus, we replace the exponential prefactor 
with general functional form, Dp G(kfi f a v ) with 
cr v being a constant, just ignoring the spatial cor- 
relation of A\. 

• In the bracket, the exponential factor e^ lAl is very 
likely to affect the structure of BAO in the power 
spectrum, since the bracket includes the term lead- 
ing to the Kaiser effect in linear regime. At the 
scale of our interest, the contribution coming from 
the factor e^ lAl should be basically small, and the 
perturbative expansion may work well. Regarding 
the quantity ji as a small expansion parameter, we 
perturbatively expand the terms in the bracket of 
the integrand. Up to the second order in ji, we 
have 

(e^A 2 A 3 ) c + (e^A 2 ) c (e^A 3 ) c 
~ (A 2 A 3 ) c +j 1 (A 1 A 2 A 3 ) c 

(A 1 A 2 ) c (A 1 A 3 ) c } + 0(f 1 ). 



+ ] 1{\{AIA 2 A 3 ) C 



(7) 



Then, from Eq. ©, the proposed model of RSD is 
expressed as 

P iS \k,fi) =D FoG [ktifa v ] 

X {P Kais er(fc,Al)+A(fc ) /i)+B(fc, ) t 1 )+£)(fc ) /i)}. (8) 

Specific form of the damping function -Df g will be given 
later [see Eq. ([34|) ]. Here, the power spectrum fkaiser is 
originated from the term (A 2 A 3 ) C , and it corresponds to 
the nonlinear generalization of the Kaiser formula fre- 
quently used in the literature (e.g., [f| EH): 

PKaiscr^, m) - Pss(k) + 2 //i 2 P S e(k) + /V Pee{k). 

(9) 

The functions P$s, Pse, and Pee are the auto-power 
spectra of density and velocity divergence, and their 
cross spectrum. The velocity divergence is defined by 
6 = —Vv/(faH). On the other hand, the terms A, B 
and D induce the corrections arising from the low-fc ex- 
pansion, and these are respectively defined by 



A{k,n) = Ji J d 3 xe ikx (Mij),, (10) 
B{k,fi) = jf j d 3 xe* kx {A 1 A 2 ) c {A 1 A a ) c . (11) 

(12) 



D(k, fi) = f / d 3 xe ik x (Ji. !>.!,) 
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The term D turns out to be of higher order if we employ 
the perturbation theory calculation, and in Ref. [24[, it 
has been neglected. While we basically follow their treat- 
ment here ignoring the D term, as already discussed in 
Ref. [24| , this is a priori no longer consistent at two-loop 
order. The validity of such an heterogeneous treatment, 
and specifically the impact of D term on the predictions 
of redshift-space power spectrum will be later discussed 
in detail (see Sec. PvT) . 

In computing each term of the expression © with 
the PT treatment of large-scale structure, we invoke a 
single-stream approximation, in which the dynamics of 
large-scale structure is described by the density S and 
velocity divergence 6. The expression in Eq. ([9]) is also 
the outcome of the single-stream approximation, and it 
is usually valid as long as we are interested in the linear 
and weakly nonlinear regimes of the gravitational clus- 
tering. Then, the expressions for the terms A and B can 
be recast as 

x {B a (p, k p, -k) - B a (p, k, k - p)} , (13) 
B(k, u) = (fe/x f f J ^F(p)F(k ~ p) ; (14) 

F(p) = ^{p se (p)+ f4Poe(p)\, 
V I P J 

where the function B a is the cross bispectra defined by 

(e{k{) |«5(fc 2 ) + / |p(fc 2 )} [s(k 3 ) + f |p(fc3)}) 

= (2^) 3 ( 5 D (fe 1 + k 2 + fc 3 ) B a {k 1 ,k 2 , fe 3 ). (15) 

Note that while we employed the low-A: expansion, we 
do not assume that the terms Ai themselves are entirely 
small. In this sense, the expression ([8]) with the cor- 
rections (fl~3|) and (fT4|) still holds some non-perturbative 
properties. A more detailed study revealed that the A 
and B terms basically give the positive contributions, 
and moderately but notably affect the shape and struc- 
ture of BAOs. In particular, as revealed by Ref. [26j], the 
A term exhibits a strong dependence on the halo/galaxy 
biasing, leading to a large-scale enhancement in ampli- 
tude relative to the real-space clustering. The effect is 
especially prominent for massive halos or highly biased 
objects, and in the presence of A and B terms, the model 
© indeed reproduces the halo redshift-space clustering 
quite well. 

In this paper, based on the T-expansion, we will make 
a fully consistent calculation of the redshift-space power 
spectrum [Eq. ©], including the PT corrections up to 
the two-loop order. With the regularization proposed in 
Ref. [i^], the power spectrum can be computed with a 
well-behaved UV behavior, which enables us to give a 
quantitative predictions for the correlation function. 

Finally, we briefly mention other improved models and 
treatments proposed recently, and note their qualitative 



differences. Ref. [28[ proposed the PT model based on 
the streaming model, allowing the scale-dependent ve- 
locity dispersion. Including the non-Gaussian correc- 
tions computed with standard PT, the model success- 
fully describes the anisotropic correlation functions. Ac- 
curate prescriptions for the anisotropic correlation func- 
tions have been also given in Refs. [13, [2{|, but these 
are constructed based on the Lagrangian PT. On the 
other hand, Ref. [3l[ presented an alternative PT pre- 
scription for redshift-space power spectrum, based on the 
moment-based expansion proposed by Ref. (30| . In this 
treatment, the higher-order corrections of RSD are all 
measurable quantities in the iV-body simulation [H, |44[ , 
and Ref. [3l[ compared their PT results with simulations 
term by term. Incorporating the effects of the small-scale 
velocity dispersions, the model has been shown to accu- 
rately describe the power spectrum in weakly nonlinear 
regime. Our RSD model may be regarded as a semi-PT 
model in the sense that a part of the terms (i.e., damp- 
ing function) is not perturbatively treated, introducing 
a free parameter. Nevertheless, with the T expansion, 
the model for the first time gives a consistent and accu- 
rate prediction for both the power spectra and correlation 
functions. 



III. THE T EXPANSION AND COMPUTATION 
OF REDSHIFT-SPACE POWER SPECTRUM 

The expression of redshift-space power spectrum in 
previous section involves not only the real-space power 
spectra but also the higher-order corrections such as bis- 
pectrum. Although our main focus is the weakly non- 
linear regime of gravitational clustering, the standard PT 
is known to produce ill-behaved higher-order corrections, 
that lack good convergence properties. Therefore, prac- 
tical calculations of redshift-space power spectrum are 
better made with a resummed PT scheme, improving the 
convergence of PT expansion on small scales, so that cor- 
relation functions can be safely computed. In this paper, 
we consider the T expansion, adopting the prescription 
for the regularization PT calculation by Ref. [4(3] . 



A. Regularized V expansion 

As seen in previous section, the density 5 and veloc- 
ity divergence 9 play an important role to describe the 
redshift-space power spectrum. Let us introduce the two- 
component multiplet: 

* a (fc) = (<5(fc), ®) , (16) 

where the subscript a — 1, 2 selects the density and 
the velocity components. For our interest of the weakly 
nonlinear scales, where the single-stream approximation 
gives a very accurate prescription, the evolution of $ a 
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is governed by the dynamics of the self-gravitating pres- 



sureless and irrotational fluid flow 41] . To perturbatively 



solve the equation for fluid dynamics, a naive treatment 
with the standard PT is to just expand the fields 4* a 
in terms of the initial fields. For the late-time epoch 
at which the growing-mode contribution is dominant, we 
then formally obtain the following expression: 

*°( fc ) = / " 3 ( t)3(nl)" Mfc - fcl-n) 

X ?(")(*!,•■■ ,fcn)«o(fcl)---«o(*! n ). ( 17 ) 

where So is the initial density field, and rj = laD(t) with 
the quantity D being the linear growth factor. The kernel 
Fa is the symmetric function, and sometimes written 
as — (F n , G n ), whose explicit expressions can be 

recursively obtained |4l| . 

As we mentioned in Sec. U however, the standard PT 
is known to produce ill-behaved high-A; behavior. This 
prevents us from obtaining a convergent result for the 
correlation function. Here, as alternative to the standard 



PT framework, we consider the T expansion, which is one 
of the non-perturbative PT frameworks. In this scheme, 
the multipoint propagator constitute the building blocks, 
and all the statistical quantities can be expressed in terms 
of these propagators. Denoting the (p + l)-point propa- 
gator by I>\ we define 



1 / SP^ a (k, V ) 

p\ \6S (k 1 )---SS (k p ) 



Sn{k - ki... f 



x (27r )3( P -i) ri p) (fci,---,fc P ;7?)- (18) 

With these objects, the power spectra are shown to be 
expressed as |38j |. 

x rW(g Xl ■ • • , g^rf • • • , q t -n) P (qi) ■ ■ ■ P (q t ). 

(19) 

Further, the bispectrum is expressed as 



B abc (k 1 ,k 2 ,k 3 ;7j) = ^ 



r + s \ ( s + t \ ( t + r 



-\s\t\ 



d 3 pi ■ ■ ■ d 3 p r d 3 q, ■ ■ ■ d 3 q„ d 3 r\ ■ ■ ■ d 3 r t r . r ,, . r . 

( 2 ^)3(r-l) ( 27r )3( S -l) ( 27r )3(«-D _ _ &(fe2 + 91 - _ ri - l) SD ^ 3 + ri ~ 4 ~ 



x T[[ +s) (p 11 ■ ■ ■ ) p r ,q 1 ,- - ,q s ;ri)T^ +t) {-q 11 ■ ■ ■ ,q s ,r ll - - , r t ; rj) r( t+r) (-ri, • • • ,-r t ,-p 1} --- ,-p r ;rj) 
x PoOl) • • • Po(Pr)Po(qi) ■ ■ ■ Po(q s )Po{ri) ■ ■ • i%(r t ). 



(20) 



The multipoint propagators are the non-perturbative 
quantities that can be obtained by summing up a class 
of infinite series of the standard PT expansion. The im- 
portant remark is that one can exploit the asymptotic 
properties of the propagator in both low- and high- 
k regimes. To be precise, in the high-A; limit, higher- 
order contributions can be systematically computed at 
all orders, and as a result of summing up all the con- 
tributions, the multi-point propagators are shown to be 
exponentially suppressed [38, 39], 



T^ikx,--- ,k p ;r,) 
— > exp 



r<W ,k P ;v) (21) 



with k = |fei H h k p \. Here, the quantity r„^ ree is the 

lowest-order non- vanishing propagator obtained from the 
standard PT calculation, and a v is the one-dimensional 
root-mean-square of the displacement field defined by, 



1 f d 3 q P (q) 
3 J (2tt)3 g 2 • 



(22) 



On the other hand, at low- A;, the propagators are ex- 
pected to approach their standard PT expressions that 
can be written formally, 



' ' ' i k p ] rf) — r^ tree (fei, • • • , k p ; rf) 

oo 

+ E r il 1 oo P (fci,---,fe P ;ry). (23) 

n=l 



For the dominant growing-mode contribution we are in- 
terested in, each correction term is expressed in terms of 
the standard PT kernels as, 



itLeCfci,"- ,k p ; V ) = e^F^(k 1 ,c...,k p ), (24) 



for the tree-level contribution, and 
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rS-ioo P (fei: ■ • • . = e {2n+p) n ^ J ^ p [ 2 ^ Pn ^i 2n+p) (Pi, -pw ,p n , -p n M, • • ■ , k p )p ( Pl ) ■ ■ ■ p Q ( Pn ) 



^a,n-loop(kli ' " ; &p) 



(25) 



r 



for the n-loop order contributions, with the coefficient 

(p) • -u 
Cn given by 



i rt = (2n-l)H. 



(26) 



(p) 

Note that n-loop order correction „_ loop is that each 
perturbative correction possesses the following asymp- 
totic form, 



r 



(p) 

a, n-loop 



l ( kW>\ r(p) 



(27) 



which consistently recovers the expression (1211) when we 
sum up all the loop contributions. This fact leads to a 
novel regularized scheme, in which the low- and high-fc 
behaviors are smoothly interpolated without any ambi- 
guities [40j |. Then, the regularized propagators are ex- 
pressed in a transparent way in terms of the standard 
PT results, and one gets 



ri P rcs-(ki, • • • , k p ; rj) — e 



v >i 



F^(k u --- ,k p ) 1 



r 



c 1 a,l-loopl K 'l 



• hp) 



exp 



(28) 



which consistently reproduces one- loop PT results at low- two- loop order PT corrections at low-fc: 
k. This construction is easily generalized to include the 



J 



ri P rcg(^lj ' ' ' > kp'i v) 



FM(k u --- ,k p )\l 



fcVV i AVV 



e 2 «r ( £ oop (k 1 ,---,kp)\l 



k 2 a 2 e 2 ^ 



5 V ^Vldc-p^li ' ' ' 7 k p ) 



exp 



r 



(29) 



— (p) 

Note that the functions r„_i oop are the scale-dependent 
part of the propagator defined by Eq. (|2"5)) . 



B. Redshift-space power spectrum from 
regularized V expansion 

In what follows, with the regularized prescription of 
the multi-point propagators, we apply the T expansion to 
calculate the redshift-space power spectrum and correla- 
tion function up to the two-loop order. We hereafter call 
the PT treatment with regularized T expansion RegPT. 
Here, we briefly describe the technical implementation of 
the RegPT to the model of RSD. Readers who are only 



interested in the results may skip this subsection and di- 
rectly go to next section. 

In Eq.© ignoring the D term, there appears three 
terms to be computed perturbatively, i.e., -Fkaiser, A and 
B terms, which include the real-space power spectrum 
and bispectrum. Below, we will separately give a pre- 
scription how to compute each term. 



1. Nonlinear Kaiser term PKaiser 

The power spectrum PKaiser consists of the three real- 
space power spectra, P55, Pse, and Pgg. In terms of the 
power spectra P a b for the doublet \l/ a , these are equiva- 
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lent to Pn, Pi 2, and P22, respectively. Thus, for a prac- 
tical computation of Pk^s C1 -, we just follow the prescrip- 
tion presented in Ref. [451 ) . and use the same technique 
to calculate each power spectrum contribution at one- 
and two- loo order. In Appendix [XJ we give explicit ex- 
pressions for the power spectra P a b up to the two-loop 
order. 



Note that the A term itself is already a higher-order 
contribution, and in computing the redshift-space power 
spectrum, the tree-level and one-loop calculations of the 
bispectrum are sufficient for a consistent calculation of 
P^ up to the one- and two-loop order, respectively. 
The expressions for the regularized bispectra are explic- 
itly given in Appendix [XJ 



2. A term 



The A term given in Eq. ([T3|) includes the cross bispec- 
trum B a [Eq. (fT5)) ]. and we thus need the RegPT to ex- 
plicitly compute this term. The function B a is expressed 
in terms of the bispectra of ^Sf a : 

B a (ki,k 2 , k 3 ) = B 2 u(kt,k 2 , k 3 ) 

+ f(j^) B 221 (ki,k 2 ,k 3 ) 

■ /(tt-) B 2 i 2 (ki,k 2 ,k 3 ) 

+ / 2 (^f ) 2 B 222 (ki,k 2 ,k 3 ). (30) 



In computing the A term, the expression given in 
Eq. (p~3|) is not suited for a practical purpose. Here, fol- 
lowing the same technique as used in Ref. [45J , we derive 
alternative expression for which the term is expanded 
as the polynomials of /i and /. The detail of deriva- 
tion is described in Appendices A and B of Ref. [2~jj . 
and we here present the final expression, in which the 
three-dimensional integral is reduced to the sum of the 
two-dimensional integrals: 



k 3 f°° f 1 r ~ 1 

Mk,ti = J2 E ^ n f a+b ~ l 7T^ / dr dx [A2 b (r,x)B 2ab (p,k- P ,-k) + A2 b (r,x)B 2ab (k-p,p,-k)j. 

n=la,b=l * 71 ' ^° 

(31) 



Here, r and x are the dimensionless variables defined by 
r = p/k, and x = (k ■ p)/{kp). The non- vanishing com- 
ponents of A bc and A bc are summarized as follows: 



3. B term 



A\ x — r x, A 21 — 



r 2 (-2 + 3rx){x 2 - 1) 
2(1 + r 2 - 2rx) ' 
_ r{2x + r(2-6x 2 )+r 2 x(-'3 + 5x 2 )} 



A21 2(l + r 2 -2rx) 

4 2 _ a\ a2 _ Al A3 _ a2 

^12 ~~ ^11' ^22 — ^21 ' ^22 — "-2\i 

2i r 2 (rx-l) ~i _ r 2 {-l + 3rx){x 2 - 1) 

n-~TTZ3 A 2 i- 2 {l + r 2 -2rx) ' 



1 + r 2 — 2rx 

r 2 {\ - 3x 2 + rx(-3 + 5x 2 )} 
2(1 + r 2 - 2rx) ' 
A 2 — A 1 A 2 — A 1 A 3 — A 2 



A 2 



The expression of the B term in Eq. (TH|) involves the 
integral of the power spectra, to which we apply the 
RegPT. Note that the B term itself is already higher- 
order contribution, and for a consistent calculation of the 
redshift-space power spectrum, the tree-level and one- 
loop calculation of the real-space power spectra are suf- 
ficient. 



In Ref. [24|, the expression of the B term suited for a 
practical calculation has been derived without employing 
the perturbative calculations: 



B(*,M) = E 



2 

E 

n— 1 a,b— 1 



M 2n (-/) 



a+b 



(2n) 2 



dr 



dxB2 b (r, x) K —— 

1 U 



2rx)P b2 (kr) 



2rx) a 



(32) 
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TABLE I: Cosmological parameters for iV-body simulations (ACDM) [45[ 



Name 




# of particles 




# of runs 






nb/fim 




ft 


n a 


CT8 


wmap5 


2,048/i _1 Mpc 


1, 024 3 


15 


60 


0.279 


0.721 


0.165 


-1 


0.701 


0.96 


0.815 9 



The coefficients B™ b are given in Appendix B of Ref. [24 



IV. RESULTS 

In this section, we present the results of the PT cal- 
culations, and compare the PT predictions with ./V-body 
simulations. The redshift-space power spectrum and cor- 
relation function are computed with the RegPT consis- 
tently including the PT corrections up to the one- and 
two-loop orders. Note that for the one- and two-loop pre- 
dictions, the bispectra and power spectra in the A and B 
term have been consistently computed including the PT 
corrections up to the tree- and one-loop orders, respec- 
tively. After briefly describing the ./V-body simulations in 
Sec. II V Ai the comparison between the predictions and 
simulations is presented for power spectrum in Sec. IIVBI 
and for correlation function in Sec. II V CI 



A. A-body simulations 

To compare the model prediction with ./V-body sim- 
ulations, we use the data set presented in our previ- 
ous paper (45|. The data were created by a public N- 
body code GADGET2 46] with cubic boxes of side length 
2,048/i~ 1 Mpc, and 1,024 3 particles. The cosmological 
parameters adopted in these ./V-body simulations are ba- 
sically the same as in the previous one, and are deter- 
mined by the five-year WMAP observations [47J (see Ta- 
ble [T| . The initial conditions were generated by the 2LPT 
code [48| with the initial redshift z; n i t = 15, and the re- 
sults of 60 independent realizations are stored at redshifts 
z = 3, 2, 1, and 0.35. The total volume of each output is 
515/i- 3 Gpc 3 . 

We measure both the matter power spectrum and cor- 
relation function in redshift space, applying the distant- 
observor approximation. For the power spectrum, we 
adopt the Cloud-in-Cells interpolation, and construct the 
Fourier transform of the density field assigned on the 
1,024 3 grids. As for the estimation of two-point corre- 
lation function, we adopt the grid-based calculation us- 
ing the Fast Fourier Transformation 35]. Similar to the 
power spectrum analysis, we first compute the square of 
the density field on each grid of Fourier space. Then, 
applying the inverse Fourier transformation, we take the 
average over realization, and finally obtain the two-point 
correlation function. The implementation of this method, 
together with the convergence test, is presented in more 
detail in Ref. 35]. In what follows, the error bars of 



the ./V-body results indicate the root-mean-square fluctu- 
ations of the averaged power spectra or correlation func- 
tions over the 60 realizations. 
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FIG. 1: Contribution of the correction terms for the PT model 
of redshift-space power spectrum. Left and right panels re- 
spectively show the A and B terms at z — 1, respectively. For 
illustrative purpose, all the results are multiplied by k' i/2 . The 
A and B terms are expanded as A[k, fi) = X]„ I- 1 - 2 ™ A2n(k) and 
B(k,[i) = 2~Dn M 2 ™ ^2n(fe). Here, we plot the scale-dependent 
coefficients A^n and i?2n (A2, -B2: magenta, A4, B4: cyan, 
Aq Bg: green, Bg: yellow). In each panel, solid lines are the 
corrections for the two-loop contributions, while the dashed 
lines are the results for the one-loop contributions. For ref- 
erences, the power spectra Pss, Pse, and Pee computed from 
RegPT are also shown in black, blue and red lines, respec- 
tively. 



B. Power spectrum 

Before presenting a quantitative comparison, we first 
look at the contribution of each term appeared in the 
model of RSD. In Fig. [T] for specific redshift z = 1, we 
plot the results for A and B terms as well as the power 
spectra Pss , Pse , and Pee ■ From the expressions given in 
Eqs. (|5H and (|3"2"j) . the A and B term can be expanded as 

#,^) = E^2„(fc)/i 2 " and B(k,fi) = E^ 2 „(fe)/i 2 «, 
and we here plot the scale-dependent coefficients A 2n and 
B 2n multiplied by fc 3 / 2 (A 2 , B 2 : magenta, A4, B4: cyan, 
A 6 , B 6 : green, B s : yellow). Dashed and solid lines re- 
spectively indicate the one- and two-loop contributions to 
the redshift-space power spectrum. For reference, we also 
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FIG. 2: Ratio of power spectra to the smoothed reference spectra in redshift space, (k)/P l ' no _ wigglc (&). Left, middle, 
and right panels respectively show the monopole (£ = 0), quadrupole (I = 2) and hexadecapole (£ = 4) contributions to the 
redshift-space power spectrum. N-body results are taken from the wmap5 simulations of Ref. [4a]. The reference spectrum 
no-wiggic ^ s computed with the no- wiggle approximation of the linear transfer function [if|, taking account of the linear 
theory of the Kaiser effect. Long-dashed and solid lines respectively indicate the results based on the RegPT calculations 
at one- and two-loop orders, adopting the Gaussian (thin) and Lorentzian (thick) form of the damping functions. Triangles 
and squares in middle and right panels are also obtained from the same calculation at one- and two-loop orders, but taking 
account of the effect of finite grid-size for the power spectrum measurement in TV-body simulations (see text and Appendix iBl in 
detail). For comparison, the 1-a statistical errors of the hypothetical survey with volumes V = 5h~ 3 Gpc 3 and number density 
n — 5 x 10~ 4 h 3 Mpc -3 are estimated from Eq. ()35[). and are depicted as green shaded regions around the N-body results at 
2 = 0.35 and 1. 



plot the results with standard PT calculations in dotted 
lines. Compared to the standard PT results, the coeffi- 
cients of the A term are slightly enhanced at the two-loop 
order, and the oscillatory feature originating from the 
BAOs is somewhat smeared. This is similar to what we 
saw in the real-space power spectrum. Fig. [T] apparently 
indicates that at the two-loop order, the A term seem to 
eventually dominate the total power spectrum at small 
scales. However, this is actually not true. Because of 
the exponential cutoff generically appeared in the multi- 
point propagators, the amplitudes of both the A and B 
terms are suppressed at small scales, as similarly seen 
in the power spectra of the density and velocity fields. 
This regularized UV property enables us to give a con- 
vergent result for the correlation function, although, as a 
trade-off, the prediction of redshift-space power spectrum 
eventually becomes inappropriate at some small scales. 

Now, let us compare the predictions of redshift-space 
power spectrum with those of the A-body simulations in 
greater detail, and investigate the extent to which the PT 
model reproduces the simulation results. Fig.[2]plots the 
ratio of power spectra to the smooth reference spectra, 
Pi S \k)/pg _ wigglc (k), where P^ wigglc (k) is the linear 
power spectrum computed with the no- wiggle formula of 
Ref. [49J. The power spectrum P} ' is the t-th order 
moment of the redshift-space power spectrum defined by 

Pf)(fe) = ^±! J dnPW(k,riV t (p), (33) 
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Z 

FIG. 3: Fitted values of a v as function redshift. Triangles and 
squares are respectively obtained by fitting the one- and two- 
loop PT predictions of the monopole and quadrupole spec- 
tra to those obtained from the A-body simulations. Open 
(filled) symbols are the fitting results adopting the Gaussian 
(Lorentzian) form of the damping function in the theoreti- 
cal predicitons [see Eq. (|34[) ]. For reference, we also plot the 
linear theory prediction in solid line. 



with Vi being the Legendre polynomials. In Fig. [21 the 
results of the monopole (£ = 0), quadrupole (£ = 2), 
and hexadecapole (I — 4) power spectra are respectively 
shown from left to right panels. Here, the dashed and 
solid lines are the results based on the RegPT calcu- 
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lations at the one- and two-loop orders. In plotting 
these predictions, the velocity dispersion c v in Eq. flSJ) is 
treated as a free parameter, and is determined by fitting 
the model prediction to the A^-body results of monopolc 
and quadrupole spectra, assuming the Gaussian (thin) 
and Lorentzian (thick) forms of the damping function 
-DfoG : 



D FoG (x) 



l/(l + x 2 /2) 2 



Gaussian 
Lorentzian 



(34) 



The predictions were fitted to the simulation results up 
to the scale indicated by the vertical arrows. Note that 
this roughly corresponds to the critical wavenumber fc cr it, 
below which the RegPT prediction in real space agrees 
with iV-body simulation at a percent- level precision [45| . 
The fitted values of a v are summarized in Fig. [3] 

Overall, the agreement between predictions and sim- 
ulations is remarkably good for the monopole and 
quadrupole spectra. In both one- and two-loop results, 
the range of agreement is almost comparable to what we 
found in the real-space comparison. This is true irre- 
spective of the choice of the damping function. In par- 
ticular, the two-loop results look very similar to what we 
obtained with the closure [3J, |35| and standard PT cal- 
culations 24|. Rather, with the RegPT, the oscillatory 



feature in A term is erased, and the contribution of the 
A term to the BAO structure is somewhat reduced. As a 
result, the total sum of each contribution closely matches 
the ./V-body results better than the previous results. 

Turning to the hexadecapole power spectra, on the 
other hand, the simulation results show somewhat noisy 
behaviors, and it is bit difficult to compare those with 
the predictions depicted as continuous lines. This noisy 
structure basically comes from the fact that the power 
spectra are measured from the grid-assigned density field. 
Because of the oscillatory feature of the Legendre polyno- 
mials Ve, the measurement of the higher- multipole spec- 
tra tends to be sensitively affected by the finite grid-size. 
In order to remedy this, we also take account of the effect 
of finite grid-size in the PT prediction, and compute the 
hexadecapole spectra in the same way as we did in the 
A'-body simulations. The detailed prescription of this 
treatment is presented in Appendix IB"! 

The predictions taking account of the finite-grid size 
effect are depicted as the triangles and squares for onc- 
and two-loop calculations, respectively. At high-z, the 
PT results faithfully reproduces the noisy behavior in 
the simulation results. In particular, at low-fc, the sym- 
bols almost overlap each other, indicating that the simu- 
lations consistently recover the linear theory prediction. 
At lower redshifts, on the other hand, there appear small 
but non-negligible discrepancies. To see the significance 
of this, we consider the hypothetical galaxy survey, and 
estimate the expected l-er statistical errors, AP^>, de- 
picted as green shaded region at z = 0.35 and 1 in Fig. [21 



(s) 

The statistical error AP f is simply computed with 

[APf(*)] a = ^(*) 5 ** = ^^S> ( 35 ) 
with the function apj given by 



,{k) 



{21 + 1) 



dfj, 



(36) 



where -P lin is the linear power spectrum, and the Ak 
is the bin width for which we adopt the same bin size 
used in the A^-body data. The volume and the num- 
ber density of the hypothetical survey are respectively 
set to V = 5/i~ 3 Gpc 3 and n = 5 x 10~ 4 /i 3 Mpc~ 3 , 
which roughly correspond to those of the Baryon Oscil- 
lation Spectroscopic Survey (BOSS) 1 or the survey with 
Subaru Measurement of Images and Redshifts (SuMIRe) 
with Prime Focus Spectrograph (PFS) 2 3 . Then, we 
found that the predicted monopole and quadrupole spec- 
tra agree with simulations well within the statistical er- 
ror, while the discrepancy in the hexadecapole spectra is 
marginal or even exceeds the l-a error at high-fc, depend- 
ing on the functional form of the damping function -DfoG • 
These results not only indicate the sensitive dependence 
of the damping function but also suggest a possible deficit 
in the RSD model when predicting the higher-multipole 
spectra. This point will be further discussed in greater 
detail in Sec. IVl 



C. Correlation function 

We next consider the correlation function. We first 
show in Fig. 0] the the large-scale behavior of the cor- 
relation functions, focusing on scales around the baryon 
acoustic peak. Left, middle and right panels respectively 
plot the results of monopole, quadrupole, and hexade- 
capole correlation functions. The predictions depicted as 
solid and dashed lines are basically obtained from the 
power spectra through the relation: 



d s) (*) 



dk k 2 



(37) 



Note again that thin and thick lines represent the pre- 
dictions adopting the Gaussian and Lorentzian damping 
functions, and we use the same velocity dispersion <r v as 
determined in the power spectrum analysis. At all red- 
shifts, the one- and two-loop results do indeed agree with 
Af-body results quite well. This is to be contrasted with 



1 http://www.sdss3.org 

2 http://sumire.ipmu.jp/en/ 

s Strictly speaking, the BOSS aims at observing galaxies at 0.2 < 
z < 0.8, while SuMIRe PFS project will observe galaxies at 0.6 < 
z < 2.4 [Hot 
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FIG. 4: Redshift-space correlation functions around the baryon acoustic peak. Monopole (£ = 0), quadrupole (i = 2), and 
hexadecapole (£ = 4) moments of the redshift-space correlation function are respectively shown in left, middle and right panels. 
Dotted lines are the linear theory predictions, while the dashed and solid lines respectively represent the results based on the 
RSD model using the RegPT up to the one- and two-loop orders, adopting the Gaussian (thin) and Lorentzian (thick) damping 
function. 




FIG. 5: Redshift-space correlation functions at small scales. Plotted results are the ratio of correlation function to the linear 
theory predictions taking account of the linear Kaiser factor, i.e., ^ S ' > ( s )/d S iin( s )- For clarity, we artificially shift the results 
at each redshift by a constant value (indicated by the horizontal dotted line). Symbols and line types are the same as those in 
Fig. [2] The green shaded regions at z — 0.35 and 1 indicate the expected 1-er error of the hypothetical galaxy survey with the 
volume V — 5 h~ 3 Gpc 3 and number density n = 5 x 1CP 4 h 3 Mpc -3 . 



previous studies neglecting A and B terms [35J , in which 
the power spectra, P$$, Pse and Pee, are computed with 
closure theory. Now with the coherent treatment with 
RegPT, the model successfully describes the correlation 
functions around the baryon acoustic peak. The results 
show that at large-scales, the choice of the damping func- 
tion hardly change the prediction, and both the one- and 
two-loop predictions almost coincide each other. 

Let us look at the small-scale behaviors beyond the 
baryon acoustic scales. Fig. [5] shows the ratio of 
the correlation functions to the linear theory predic- 
tions, ^ S ' 1 (s) / ^J^j n (s) , specifically focusing on the scales 



10 h^ 1 Mpc < s < 80 hr 1 Mpc. Note that the linear the- 

(S) 

ory prediction Q y! is made with the linear power spec- 
trum taking only account of the linear Kaiser effect. As 
references, we also consider the hypothetical galaxy sur- 
vey, and show the 1-a statistical errors at z — 0.35 and 1, 
depicted as green shaded region. This is estimated from 

[Ad S) (*)] 2 = f / ^ln(ks)* P Ak)} 2 (38) 

with op } t denned in Eq. (l36|) . Here, we adopt the same 
survey parameters as we considered in Fig. [5J 

As anticipated from the power spectrum results, the 
predictions for both the monopole and quadrupole mo- 
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ments perfectly agree well with simulations at the scales, 
s > 10 — 30/i _1 Mpc, depending on the redshift. The 
range of agreement with TV-body simulations is compara- 
ble to the one obtained in real space, and roughly matches 
the range inferred from the power spectrum results. One 
noticeable point in the prediction of correlation function 
is that even the one-loop results do give an accurate pre- 
diction over a wide range of correlation function, where 
the choice of damping function hardly change the results. 

On the other hand, similar to the power spectrum anal- 
ysis, the measured hexadecapole moment suffers from the 
effect of finite grid-size, and in order to make a fair com- 
parison, we need to incorporate the effect of this into the 
theoretical calculation. The triangles and squares are 
the results taking account of the finite grid-size based on 
the prescription in Appendix [Bj Then, the predictions 
at high-z reproduce the ./V-body results almost perfectly, 
while we find a systematic discrepancy at low-z, where 
the results also show a sensitive dependence on the choice 
of the damping function. However, we note that the dis- 
crepancy seen in the correlation function is smaller than 
the statistical errors of the hypothetical surveys, and it 
seems less significant compared to the results in power 
spectrum. This is partly because many Fourier modes 
can contribute to the correlation function, and they help 
to mitigate the significance of the discrepancy seen in 
the power spectrum at some specific modes. Strictly 
speaking, the amplitudes of the correlation function are 
strongly correlated between different scales, and the error 
covariance of the correlation function may be important 
to rigorously judge the significance of the discrepancy. 
Rather, a conservative and generic statement is that for 
the scales where both the one- and two-loop predictions 
coincide each other, the non-linear effect of RSD and 
gravity can be small, and thus the discrepancy between 
the PT predictions and simulations is insignificant. 
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FIG. 6: Same as in Fig.[TJ but the scale-dependent coefficient 
of the D term at z = 1 is plotted (D2: magenta, D4: cyan, 
Dq: green, Dg: yellow). For reference, the power spectra Pss, 
Pse, and Pee computed from the RegPT at two- loop order are 
also shown. 



as 



d 3 pd 3 q p z q z 



x {Ta(p, q,-k-p-q,k)- T a {p, q,k-p,-k- q)} , 

(39) 

where the function T a is the cross trispectrum defined by 



h 2 W P 

3z afu \ I J xru ^ 1 t 4z , 



0(fci)(0(fc 2 ) \ S(k 3 ) + f -§tO(k 3 ) i i 6(k 4 ) + / -jf 0(fe 4 ) 



(27r) 3 <5 D (fci + fc 2 + k 3 + fe 4 ) T a (kx,k 2 , fe 3 , ^4)- 



(40) 



V. IMPACT OF HIGHER-ORDER 
CORRECTIONS 



So far, we have neglected the contribution of the D 
term [Eq. (|12p] to the predictions of rcdshift-space power 
spectrum. Strictly speaking, this is no longer valid when 
we perform the two-loop calculations. The higher-order 
corrections of the RSD such as the D term may help us to 
explain the small discrepancy seen in the hexadecapole 
power spectrum. In this section, we quantitatively com- 
pute the D term, and discuss the validity of the anal- 
ysis in previous section. Further, using the measured 
power spectrum of the iV-body simulation, we also es- 
timate the magnitude of a possible systematics coming 
from the higher-order corrections to the RSD. 



In deriving the above expression, we have used 
the symmetric properties of the trispectrum T a , 
i.e., T a (ki,k 2 ,k 3 ,k i ) = T a {—k\, -k 2 , -k 3 , = 
T a (fc2, fci , fc3, k*) = T CT (fci, fe2, A:4, k 3 ). The trispectrum 
T a is rewritten with the sum of the cross power spectrum 
T a bcd for the doublet ^ a as 

T a (ki,k2,k 3 ,ki) = T22ii(ki,k 2 ,k 3 ,k 4: ) 

2 



■!-/ ( ^ ) T 22 - Il {k 1 .k 2 .k,;.k i ) 



H- ./' ( ^ ] T^A^.ki.kn.kO 



k 3 k± 



T 2 222(k 1 ,k 2 ,k 3 ,k 4 ), 
(41) 



Let us first consider the D term. Eq. (|T2|) can be recast 



which can be computed with RegPT. For the predictions 
of redshift-space power spectrum at two-loop order, it is 
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FIG. 7: Ratio of power spectra to the smoothed reference 
spectra in redshift space, P£ S \k)/P^ o _ wigslc (k) for the pre- 
dictions including the D terms (blue dot-dashed). The results 
at z = 1 are specifically shown. For comparison, the two-loop 
results ignoring the D terms are also shown (magenta solid). 
Thin and thick lines are the results adopting the Gaussian 
and Lorentzian damping function, respectively. 



sufficient to give the tree-level results for T abcd , whose 
explicit expression is given in Appendix IA 31 

Fig. |6] shows the D term computed at specific redshift 
z = 1. Similar to the A and B terms, the D term can 
be expanded in powers of fi [see Eqs. (|39t and (|4Tj) ]: 
D(k,fi) = Yl=i M 2 ™ D 2n (k). We then plot the coeffi- 
cients I?2n as function of wavenumber. Note that unlike 
the A and B terms, we were not able to derive a simpler 
expression like Eqs. (|3"Tj) and (|3"2j). Hence, we employed 
the Monte Carlo technique to directly perform the six- 
dimensional integral in Eqs. ([39]) . and obtained the result 
in (fc, /i)-plane. For each wavenumber, we applied the 
multipole expansion, and characterize the /x-dependence 
of the D term by the Legendre polynomials. Finally, 
the resultant coefficients are translated into those in the 
power- law expansion, Z?2n- 

The D term shown in Fig.[6]has the amplitudes roughly 
comparable to the B term. However, a careful look at the 
/i-dependence reveals that the coefficients D2 and -D4 are 
rather small. Also, while the amplitude is non-negligible 
for Dq and D^, their signs are opposite each other. These 
facts imply that when convolving with the damping func- 
tion -DfoG j the contribution of the D term becomes neg- 
ligible for the monopole and quadrupole spectra, and is 
largely canceled for the hexadecapole power spectrum. 
This is indeed manifested in Fig. [7j where we compare 
the model prediction including the D term (blue dot- 
dashed) with TV-body simulations. Except for a slight 
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FIG. 8: The ratio of power spectra measured in (fc, fi) space 
from the simulations to those from the theoretical calculation 
[sec Eq. (|42[l ]. The results at z = 1 are plotted as function of 
kfi. The color scale indicate the power spectrum contribution 
at different range of wavenumbers. Left panel shows the re- 
sults ignoring all the corrections, A, B and D terms, in the 
theoretical calculation. The middle and right panels are the 
results taking account of the corrections, but in middle panel, 
the D term is neglected. For reference, we also show the 
Gaussian (£>f g = exp{— (k \i f ov) 2 }, solid) and Lorentzian 
(7?FoG = 1/{1 + (k [i f ov) 2 /2} 2 , dashed) damping functions 
adopting the fitted values of a v . 



change in the hexadecapole, which makes the prediction 
slightly better if we adopt the Lorentzian damping func- 
tion, the resultant power spectra are hardly affected by 
the D term. The fitted value of the parameter cr v almost 
remains the same. Therefore, we conclude that the ac- 
tual contribution of the D term is less significant for the 
prediction of redshift-space power spectrum. 

Nevertheless, this does not prove that the corrections 
higher than the D term originated from the low-fc ex- 
pansion in Eq. (JJJ are entirely negligible. Further, as 
we saw in previous section, there appears a small but 
non-negligible difference in the hexadecapole spectrum 
between the model predictions and TV-body simulations 
(Fig. [5]). Since a part of our RSD model has not been 
derived by the perturbative expansion, there would exist 
some missing corrections that can systematically affect 
the prediction of higher-multipole power spectra. 

In order to elucidate the potential systematics, we mea- 
sure the power spectrum of TV-body simulations in (/c, //) 
space, and divide it by the PT prediction in the following 
way: 



o(S) 
N-body 



(fc,/i) 



Pi 



Kaiser 



(fc, fi) + A(k, fi) + B(k, ft) + D(k, fi) 



(42) 



In Fig. HI the results at z = 1 are plotted as the func- 
tion of the single variable, kfi. The color scales repre- 
sent the contributions from the different wavenumbers. 
The plotted results are the residual contributions of the 
RSD that is not described by the low-fc expansion, but 
is rather characterized by the damping function -DfoG m 
our PT model [Eq. ©]. Compared to the case ignor- 
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ing the correction terms (left panel) , the residuals shown 
in the middle and right panels are mostly character- 
ized by the single-valued function of k/j,, and within the 
plotted range, they are approximately described by the 
Gaussian (JDfoG = e~ ( - fcA1 ^' Tv - ) , solid lines) or Lorentzian 
(Dp Q = 1/{1 + (fc/i/er v ) 2 /2} 2 , dashed lines) damping 
function. This is indeed what we expected from our RSD 
model. 

However, a closer look at the scatter around the damp- 
ing functions reveals some systematics that the \ow-k 
residuals tend to take larger values, while the high-fc 
residuals take slightly smaller values than the Gaussian 
or Lorentzian damping function, indicating the imper- 
fect modeling of RSD. These systematics in the scatter 
would be certainly the source of the discrepancy seen in 
the hexadecapole power spectrum. In this respect, a bet- 
ter description for the power spectrum suppression might 
be crucial for an accurate prediction of higher-multipole 
spectra. Recent study by Ref. [31( suggests that a simple 
prescription with the function -DfoG might be a oversim- 
plified treatment, and the suppression of power spectrum 
cannot be controlled by a single parameter er v . These 
two facts may pose a caution for the assumption of our 
model, and the suppression and modulation of the power 
spectrum might not be treated separately. Thus, in gen- 
eral, a more elaborate modeling of RSD is needed for an 
accurate prediction of higher multipole power spectra, 
although the actual impact of such a modeling on the 
cosmological application is unclear, and is beyond the 
scope of this paper. 



VI. CONCLUSION 

In this paper, we have implemented the resummed PT 
scheme called RegPT, which based on a multi-point prop- 
agator expansion, together with an improved model of 
rcdshift-space distortions (RSD), in order to compute the 
matter two-point statistical properties in redshift space, 
both in Fourier and in configuration spaces. For the 
first time, calculations consistently include PT correc- 
tions up to the two-loop order. Our analytical predic- 
tions are then compared with iV-body simulations. We 
found that, whenever predictions are at percent level in 
real space, a similar performance can be obtained in red- 
shift space. With a full implementation of the RegPT 
scheme, the present work puts forward proper PT cal- 
culations in redshift space, allowing the predictions to 
be greatly improved in the sense that we are now able 
to give a consistent prediction not only for the power 
spectrum but also for the correlation function. This is 
in marked contrast to the previous analysis partly using 
the standard PT treatment. 

We have also investigated the influence of higher-order 
corrections in redshift space. While our RSD model has 
been derived based on the low-fc expansion, the previous 
studies have missed the non-Gaussian correction at the 
trispectrum order called D term, whose contribution is 



comparable to the two-loop order. We have computed 
the D term explicitly, and confirmed that the amplitude 
of the coefficients given in powers of \x are comparable 
to that of the B term. The actual contribution to the 
power spectrum, however, turns out to be small and it 
would only affect the higher-multipole spectra. Never- 
theless, the present RSD model seems to have a slight 
deficit in the prediction of the higher multipoles, and 
we have quantified the validity of our model assumption 
with a help of iV-body simulations. The results in Fig. [5] 
indicate an improper treatment for the damping effect 
of the power spectrum, and this would lead to the small 
discrepancy of the higher multipoles between predictions 
and simulations. 

Nevertheless, the actual impact of the improper mod- 
eling of RSD might be less significant, considering the 
other big systematics such as the galaxy bias. Indeed, 
one crucial aspect of the galaxy bias in redshift space 
is that it can affect both the density and velocity fields, 
and the clustering feature in redshift space cannot be 
straightforwardly understood from the real-space cluster- 
ing. We have previously shown that our model of RSD, 
combined with the simple linear bias, can successfully de- 
scribe the redshift-space halo clustering in iV-body simu- 
lations quite well. Still, however, this should be regarded 
as an idealistic example. In reality, galaxies do not faith- 
fully trace the halo clustering, and peculiar velocity of 
galaxies certainly differs from the center-of-mass vel ocity 
of halos. As it has been recently advocated in Ref. |5l| . 
off-centered galaxies can have a large peculiar velocity 
due to the virial motion, and the presence of these galax- 
ies can drastically change the power spectrum (52|, [53| . 
The significance of this effect may sensitively depend on 
the type of the targeted galaxies. A careful study of the 
galaxy samples is quite crucial for the cosmological anal- 
ysis to get an unbiased cosmological constraint. 
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Appendix A: Explicit expressions for regularized 
power spectrum, bispectrum and trispectrum 



1. Power spectrum 



In this Appenix, based on the regularized T expansion, 
we give full expressions for the real-space power spectra 
and bispectra used to compute the redshift-space power 
spectrum and correlation function. 



According to the prescription described in Sec. IIII1 the 
power spectrum at one-loop order of the T-expansion be- 
comes 



/t3 
^3 lt% g (q, k g; rfiT^q, k - q; V )P (q)P (\k - q\) 



(Al) 



with the regularized propagators Treg and Treg valid at one-loop order being: 



r£e g (fc;7?) = e" 



fcV 2 e 2 " 



1 + ^_ + e2 , r (i)_ ioop(fc) 



exp 



k 2 a 2 e 2r > 



r« eg (g, k - q; r,) = e 2 ^ 2 > k - q) cxp {-^f^} 



(A2) 
(A3) 



r 



On the other hand, if we include the next-to-leading tant expression becomes 
order contribution, i.e., two- loop corrections, the resul- 



P ab (k; n) = rW cg (fc; V )r { b % g (k; V )P (k) + 2 J ^ r( 2 ) cg (g, k-q; rj)Tf} es {q, k-q; r,)P (q)P (\k - q 



d 6 pd 3 q 



ri 3 Lg(p, q,k p q- r?)rg eg (p, q,k-p-q; 7 1 )P ( P )P (q)P (\k - p - q\) (A4) 



with the regularized propagators given by 



ri reg(^! v) — e 1 



r i 2 reg(?> k - '7) = 



k 2 (T 2 e 2v 1 fk 2 a 2 e 2r > 



e 2j? ri;|_ loop (fc)|i + 



k 2 a 2 e 2r > 



e4v r a! 2-ioo P (^) 



exp 



,2>i 



, 2U , , L fc 2 a 2 e 2i n 2r ,-(2) . , .1 f fc 2 a 2 e 2 "l 
^i 2) (<7,fc-g)|l + ^— | +e 2 Tlj_ loop (g,fc-g)j cxp j 2 j 



r^L(P, 9, * - P -«;»?) = e 3 " ^J 3) (p,q,k p q) cxp J - 



fcV 2 e 2 " 



fc 2 f 7 2 e 2 " 
2 

(A5) 
(A6) 

(A7) 



where the quantity r^ oop is defined in Eq. flU . The the integrals are performed with the method of Gaussian 



higher-order contributions up to the two- and one-loop 
order of the propagators are respectively included in 
the expression of the regularized propagators ri*reg and 

(2) 

Ia,rcg, consistently with the T-expansion at two-loop or- 
der. The detailed prescription on how to efficiently com- 
pute each contribution is described in Ref. (45J. Shortly, 
we use the analytic expression of the kernels to com- 
pute r£?i oop (Ti 2 j oop for [Ij and r^ oop f or [54]), and 



quadrature. For the integral in the last term of Eq. (|A1[) . 
we adopt the Monte Carlo technique using the publicly 
available library called CUBA [55| . 



2. Bispectrum 

The real-space bispectra are needed to compute the 
A term, and for a consistent calculation of the redshift- 
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space power spectrum up to two-loop order, we need the 
bispectra at one-loop order. In this case, we can decom- 
posed the bispectrum into three pieces: 

B abc (k 1 ,k 2 , fc 3 ) = £ai, c (fci, k 2 , fe 3 ) 

+ BiHkr, fc 2 , fc 3 ) + BilUkx, fc 2 , k s ). (A8) 



These three contributions are expressed in terms of the 
multi-point propagator 10 : 



^ 6c (fe 1 ,fe 2) fe 3 ) = 2 {ri% g (fe 2) fe 3 )rg eg (A ;2 )r« g (A ;3 )Po(fe)Po(fc3) 

+ri 2 r ) cg (fei, fc 3 )r<;] eg (fci)r« g (fc 3 ) PoCWCfca) + r$ eg (fci, fe 2 )r» g (fci)r£ g (fe) Po(fci)Po(fc 2 )} , (A9) 

B&(fcl, fe 2 , fca) = 8 y it^fcx - q, q) rg og (fc 2 + q, -q) T% cg (-q - k 2 , -ki + q) 

x Podfei - q|)P (|fc 2 + q\)Po(q), (A10) 
Blll{k x ,k 2 , fca) = 6 y ^3 { r ilg(" fc 3, -fca + 9, -g) rg eg (fe 2 - g, g) r« eg (fc 3 ) P (|fe 2 - q|)P ( g )P (fe 3 ) 

+ ri 3 ) cg (-fc 2 , -fc 3 + q, -q) r« g (fe) r$ e g (fc 3 - q, q) P (|fc 3 - q|)Po(g)Po(A; 2 ) 
+ ri% g (fci - q, q) rg eg (-q, + q, -fc 3 ) r« g (fc 3 ) P (|fci - q\)P (q)Po(k 3 ) 

+ r^Lg(fcx) ifre g (-fci> - fe 3 + g, -g) r^ eg (fc 3 - g, g) P (fci)Po(g)Po(|fc 3 - 



(|fcs-g| 

q|)Po( g )P (fc 2 ) 



- rl 2 L g (fci - q, g) rg cg (/c 2 ) r^) eg (_ fel + q , q , - k2 ) p (\ kl 

+r« eg (fci) rg eg (fc 2 - q, q) r^ eg (-fc!, -q, -fca + q) Po(fci)Po(q)Po(|fc 2 - q|). } 



(All) 



For the bispectrum B 1 ,, the regularized multipoint 
propagators valid at one-loop order are computed with 
Eq. for ri^og and Eq. for r£% g . On the other 

hand, for the contributions P^L ano - ^all^ the tree- level 
propagator with exponential cutoff is sufficient for a con- 
sistent calculation at one-loop order. Explicitly, it is 
given by 

itUfci, ■■■,K) = Pi")(fc l5 • • • , k n ) exp 

(A12) 

with fc = |fci + • • • + k n \. 

Finally, for the redshift-space power spectrum at one- 
loop order, the relevant contribution in the A term is the 
tree-level results of the bispectra. The tree-level bispec- 
trum is computed from the first term in Eq. (|A8|) ■ with 
the regularized propagators given by Eq. (1A12[) . 



3. Trispectrum 

Since the lowest-order contribution to the trispectrum 
is already comparable to the two-loop corrections for the 
redshift-space power spectrum, it is sufficient for our case 
to derive the tree-level expression for the trispectrum. 
The lowest-order expression for the trispectrum T abcd be- 
comes 



T abcd (k 1 ,k 2 , fe 3 , fc 4 ) = Tl bcd (ki,k2, fc 3 , fe 4 ) 

+ T I ab I cd {k u k 2 MM) (A13) 



with the contributions T^ bcd and T^ bcd respectively given 
by " 
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T abcd( k l, k 2, k 3i k i) 



= 4 [r^) eg (fe 13 


-fe3)rg eg 


(~fel3 


-fc 4 )fl)(fel3)fl)(fc8)il(fc4) 


+ ri 2 reg(^14, 


— fc4)r breg 


(-fel4 


-k 3 )P {k 14 )P (k 3 )P {k 4 ) 


+ r^ eg (fei 2( - 


-fc 2 )r( 2 r ) cg (- 


-^12, - 


-fc 4 )Po(fcl2)Po(fc 2 )Po(fc4)- 


~ ri 2 re g (fcl4, " 


-fc4)r< 2 r ) cg (- 


-fei4, - 


-k 2 )P (k u )P (k 2 )P (k 4 ) 


+ r$ eg (fc 13> - 


-fe3)rg eg ( 


-&13," 


-fe)ft(*13)fl)(*3)fl)(*2)- 


Hri% g (fc 12 ,- 


-fc 2 )rg eg ( 


-fcl2, - 


-k 3 )Po(k 12 )Po(k 3 )P (k 2 ) 


+ r( 2 ) cg (fc 13 ,- 


"^l)r'6,reg(" 


-fel3," 


-fc 4 )P (fci 3 )Fo(fci)Po(fc4)4 


-r( 2 r ) eg (fe 34! - 


-fc 4 )r[ r e g (- 


-^34," 


-fei)Fo(A: 3 4)Po(fci) J Po(fc4) 


+ r d 2 rog( fc 34,- 


-fe 3 )ri 2 r ) cg (- 


-fe34," 


-fel)Po(fc34)fo(fc3)fo(fcl) ^ 


-rg eg (fe 14) - 


-fci)r^ r g g (- 




-fc 3 )F (fci 4 )Po(fc 3 )Fo(fci) 


+ r( 2 ) cg (fc 13 ,- 


-fei)rg eg (- 


-fel3,- 


-fc 2 )i > o(*l3)i , o(fci)i , o(*a) " 


-r( 2 r ) eg (fc 23j - 


-fc 2 )rg eg (- 


-^23," 


-^1)^23)^1)^2) 



(A14) 



and 



T a Y cd (fc 1 ,fe 2 ,fc 3 ,fe 4 ) = 



r^) eg (fe 2 , fe 3 , fc 4 ) J Po(fc 2 ) J Po(fc3)-Po(fc4) + rg^fci, fe 3 , fc 4 ) J Po(fcl)Po(fc 3 )-Po(fc4) 



(A15) 



(n) 

with the tree-level multipoint propagator T a, rcg given m 
Eq. (IA12|) . Note that the expressions for the trispec- 
trum given above reduces to those of the standard PT at 
tree order if we just neglect the exponential cutoff in the 
multi-point propagators (e.g., [56j]). 



Appendix B: Effect of finite grid-size in measuring 
the power spectrum and correlation function 

Here, we present the prescription on how to incorpo- 
rate the effect of the finite grid-size into the theoretical 
prediction of power spectrum and correlation function. 

In TV-body simulations, the measurement of the multi- 
pole power spectra is done with the density fields assigned 
on grids in Fourier space. On each grid, we first evaluate 
the square of density field multiplied by the Legendre 
polynomial, and take an average over the grids within 
the thin spherical shell around fc, the width of which is 
given by Afc. Thus, for a fixed Afc, the number of grids 
falling into the spherical shell inevitably decreases with 
decreasing wavenumber, finally leading to an inhomoge- 
neous sampling. To mimic this effect in the theoretical 
predictions, we prepare the same grid space as done in 
the analysis of the iV-body data, and assign the theo- 
retical power spectrum, given as function of fc|| and fc|| , 
on these grids. Then, multiplying Legendre polynomials, 
we take an average over the spherical thin shell. This is 
expressed as 

2£+l 



grid size is chosen as 27r/Lbox with box size Lbox = 
2, 048 hr 1 Mpc, and the width of the Fourier bin in mea- 
suring power spectra are set to Ak = 0.005 /iMpc -1 . 
Note that while the average over the spherical thin shell is 
taken over the domain with positive wavenumber, k > 0, 
we must take care when we sum up the contribution from 
the grids on the boundary of the quadrant. To avoid the 
double counting, the weight factor should be appropri- 
ately multiplied. For instance, the factor 1/2 (1/4) is 
multiplied when one (two) of the components in wave 
vector vanishes. 

In similar way, the effect of finite grid-size for the cor- 
relation function can be incorporated into the theoretical 
prediction. Note that this effect is prominent only when 
we adopt the grid-based calculation of the correlation 
function. We express the multipole moments of the cor- 
relation function with the discrete sum over the grids on 
the configuration space: 



2£ + l 



E 



S|e[si-As/2,Si+As/2] 



f W (S||,*L)7MV a )» 

(B2) 



where the redshift-space correlation function £' s ) is cal- 
culated with 



€ (S) (sii,^) 



E 

\k\e[ki-Ak/2,ki 



P (S) (fc||,fc ± )P £ (fc||/fc), 



d 3 k 

dk\\ 
2^2 



-Afc/2] 



e lk - s P^(k h k±) 

dk ± k ± P (s \k h k ± ) 

x cos(fc||S||) Jo(k±s±). (B3) 



(Bl) 

where the quantity N/. is the number of grids falling 
into the spherical thin shell. In the present paper, the 



For the results shown in Fig. [5J both the grid size and 
the width of bins As are set to 2 hr 1 Mpc in configuration 
space. 
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